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Abstract 

We consider the motion of a Brownian particle moving in a potential field and driven by di- 
chotomous noise with exponential correlation. Traditionally, the analytic as well as the numerical 
treatments of the problem, in general, rely on Fokker-Planck description. We present a method 
for direct numerical simulation of dichotomous noise to solve the Langevin equation. The method 
is applied to calculate nonequilibrium fluctuation induced current in a symmetric periodic poten- 
tial using asymmetric dichotomous noise and compared to Fokker-Planck-Master equation based 
algorithm for a range of parameter values. Our second application concerns the study of resonant 
activation over a fluctuating barrier. 
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I. INTRODUCTION 



Langevin equationi^ is the key stochastic differential equation that encompasses a wide 
range of physical, chemical and biological sciences. Although originally envisaged as a de- 
scriptor of the motion of a Brownian particle in a fluid in thermal equilibrium, the equation 
is now the basic paradigm for almost any physically realizable stochastic process which 
takes into account of the motion of particle in phase space and at the same time acted upon 
by a noise force which may be internal or external depending on whether the fluctuation- 
dissipation relation is satisfied or not, respectively. Over the years various noise processes 
with arbitrary correlation time have been in use in a wide variety of problems, e.g., in ac- 
tivated processes in kinetics 4 ^, noise-induced transport processes ^ 7 ' 8 ! 9 ! 10 ! 11 ! 12 ! 13 ! 14 ! 15 ! 16 ' 17 in 
condensed matter and biological physic o 18 ^ 19i2( | 1121 etc. Among others dichotomous or tele- 
graphic nois o 22 i 23 i 24 i 25 i 26 i 27 i 28 plays the typical theoretician's paradigm of fluctuations in this 
context. While for linear potentials the Langevin equation with dichotomous noise is solv- 
able analytically, arbitrary nonlinear potential poses serious problems. For example, it is 
very difficult, if not impossible, to treat a ratchet device with periodic potential and finite 
correlation time of noise. One therefore routinely takes resort to an equivalent Fokker- 
Planck approach and the algorithm based on this description. The notable point is that the 
approach bypasses the generation of dichotomous noise as such and a close survey of the 
recent literature suggests that a direct simulation of Langevin dynamics of a particle driven 
by dichotomous noise still remains to be addressed. Our object in this paper is thus twofold, 

(i) We present a general numerical algorithm for generation of dichotomous noise with 
exponential correlation and numerically simulate the Langevin dynamics for the particle 
without taking resort to any Fokker-Planck or master equation description. 

(ii) As two prototypical applications with additive and multiplicative noise processes we 
consider dichotomous noise-induced transport or current in a periodic potential in a ratchet 
device and resonant activation over a fluctuating barrier. 

In implementing the scheme it is essential to take care of the fact that the numerical error 
should not bring any additional tilt to the potential or break of any inversion symmetry 
or detailed balance of the system. Second, the forcing terms must be unbiased so that 
after appropriate averaging over time, space and ensemble no directed transport should 
appear as an artifact. These considerations are particularly important for thermodynamic 
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consistenc y 29 i 30 i 31 i 32 i 33 i of the numerical scheme. Thus although the numerical simulation 
of Langevin dynamics with dichotomous noise is useful in a wide variety of other situations, 
our choice of application here is guided by these considerations. 

The outline of the paper is as follows: In Sec. II we present our method of numerical 
solution of Langevin equation in presence of an internal noise and a dichotomous noise. 
Although the elements of the present numerical scheme are commonly used in Monte Carlo 
calculations, the implementation in generating dichotomous noise in the numerical solution of 
Langevin equation for realization of stochastic path and nonequilibrium fluctuation-induced 
transport and resonant activation are new in the current context. Sec. Ill is devoted to the 
numerical results on the noise induced current for a cosine potential. We compare the meth- 
ods based on Langevin equation and Fokker-Planck equation^ 3 - for a range of parameters. 
Another application is the study of activated escape over a fluctuating barrier. The paper 
is concluded in Sec. IV. 

II. NUMERICAL SIMULATION OF LANGEVIN DYNAMICS WITH DICHOTO- 
MOUS NOISE 

In order to motivate our numerical simulation method we consider the motion of an 
overdamped particle in a periodic potential, V(x), simultaneously subjected to an internal 
thermal noise, £(£), and an external dichotomous or telegraphic noise, rj(t). This is described 
by the Langevin equation of the following form 

x(t) = -V'(x)+C(t)+v(t) (2.1) 

The periodic potential V(x) may be symmetric or asymmetric depending on the specificity 
of the situation. £(£) is the thermal, Gaussian, white noise whose mean and variance are 
given by 

m) = (2.2) 
mm) = 2k B T5(t-t') (2.3) 

respectively, where ks is the Boltzmann constant and T is the absolute temperature 
which is the measure of the strength of internal noise. 
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rj(t) is the dichotomous noise which can assume only two random values, say, a and b. 
We also require the random number sequence to satisfy: 

(V(t)) = (2.4) 

<i7(*M0> = 7 gx p(-^t^) ^ 2 - 5 ) 

or in other words , the dichotomous noise must have a zero mean and be exponentially 
correlated. Q is the strength and r is the correlation time of the dichotomous noise. 

Our problem here is to find out the sample path x(t) of the particle as a function of time 
and calculate the relevant quantities for nonequilibrium fluctuation- induced transport. More 
specifically, we calculate the average velocity of the particle under a steady state condition. 

A. The algorithm 

A simple approach for numerical simulation of Eq. (|2.1jl is to discretize time t and to use 
a predictor-corrector method in advancing the particle from x(t n ) to x(t n+ \) as follows: 

xi(t n+1 ) = x(t n ) - V'{x{t n ))At + rj(t n )At + {2k B TAt) 1 ' 2 W n (2.6) 

x(t n+1 ) = x(t n ) - ~ [V'(x(t n )) + V^xiftn+i))] At + v (t n )At + {2k B TAt) 1 ' 2 W n (2.7) 

where the first step ()2.6|) implies the simplest Euler-type algorithm for the predictor, 
whereas the second step corresponds to the corrector. At is the time step and W n is a Gaus- 
sian distributed random number with zero average and a variance of unity independently 
chosen at each step using a standard Box-Muller algorith m 35 i 36 i 37 i 38 . 

The next important step is the generation of dichotomous random number rj(t) with 
given properties (|2.4l 12. 5J) such as, zero average and exponential correlation. To this end we 
proceed as follows: 

B. Generation of dichotomous noise 

We consider a random variable rj(t) which switches between two values a and b randomly 
in time. The rate of switching from a to b is fj, a and from b to a is This two-step process 
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can be described by a probability loss-gain equation or master equation 1 



—P(a,t\x,t ) = -fx a P(a,t\x,t ) + fi b P(b,t\x,t ) (2.8) 
at 

j t P{b,t\x,t ) = fi a P(a,t\x,t ) - n b P(b,t\x,t ) (2.9) 

where P(a,t\x,to) is the conditional probability that the variable rj(t) will assume the 
value a at some time t given that it was x at earlier time t . P(b,t\x,t ) can be defined 
similarly. 

The conservation of total probability demands 

P(a,t\x,t ) + P(b,t\x,t ) = 1 (2.10) 
The initial condition for Eqs. (j2.8|) and ()2.9j) are given by: 

P(x',t\x,t ) = 5 x , x at t = t (2.11) 
Eq. ()2.8|) and ()2.9|) admit the following steady state solutions 



P s (a) = P(a, oo\x, t ) = - (2.12) 
P s (b) = P{b,oo\x,t ) ''" 



The conditional probabilities P(a,t\x,t ) and P(b, t\x, to) for time t can be obtained by 
solving the Eas.(l2~£l) and flXHJ) subject to the conditions (|2~Pnj) . (EHTl) and and (l2~T2l : 



P(a,t\x,t ) = — — h ( — — — Sax — — S bx ) exp [-(/i a + /i fe )(t - t )] (2.13) 

Va + Hb KHa + ^b Ma + A*6 



P(MMo) = — v— - ( — v— ^ - — — J exp [-(/i a + ^)(t - t )] (2.14) 

Ha+^b \t L a + fJ-b /Ja+yUfe / 

Eqs. (|2.13|) and ()2.14|) can be used to calculate the mean and variance of the dichotomous 
noise in the steady state as 

(rj(t)) = afib - bfia (2.15) 
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m {()) = f m^.\' + ^(a + »)' ^ + _ 

respectively. Eqs. (j2.15j) and (J2.16j) can be identified as dichotomous noise with given 
properties (|2.4|) and (|2.5jl . respectively, provided 



with the following identification 



afi b -bfi a = (2.17) 



t= (2.18) 

/J-a + fJ>b 



and 



Q fJLafJLb{a + bf 



(2.19) 

Now Eq. (j2.17|) can be rearranged to show — a ^ which when used in fl2.19|) 



7" (A*a + ^b) 2 

{fla+IJ-b) 2 



yields 



Q = abr (2.20) 
Furthermore, for convenience, we define an asymmetric parameter 

0=|a|-|6| (2.21) 

For choosing a set of seven parameters a, b, /i a , Q, 9 and r it is pertinent to satisfy 
four relations (|2.17|) . (|2.18j) . ()2.20|) and ()2.21|) . This implies that only three independent 
parameters are needed to specify the dichotomous noise rj (t) with zero mean and exponential 
correlation. 

With these aforesaid preliminaries we are now in a position to generate realizations of 
a stochastic process for a dichotomous noise. This is done in the following way. Let the 
particle is located initially (t) at x n = a. To determine whether the particle moves at time 
t\ = t + At to another site x n +i = b or remain at the same site x n = a, we consider the 
following conditional probability as given by (|2.1H|) 

P(a,t 1 \a,t) = ^ + ^ exp[-(/i a + ^)At] (2.22) 

A*a + ^b L*>a + H 
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An uniformly distributed random number R between [0, 1] is now generated by the 
computer. This number is then compared against the conditional probability (|2.22|) . If 
P(a,ti\a,t) > R, then we accept the value of the noise a, i.e., x n = a, else we accept the 
value b or x n+ \ = b. If the value of noise is b at £i(= t + At) then we calculate the conditional 
probability of jumping to another site x n+2 = a at t 2 (= t\ + At) as 

P(a,t 2 |Mi) = — r ^— exp[-(^ + ^)At] (2.23) 

Ma + lib Ha + Hb 

On the other hand if the value of the noise is x n+ i = a at £i(= t+At) we calculate the con- 
ditional probability P(a, tzla, t\) using ()2.22j) . We then compare the probability P(a, t 2 \b, ti) 
or P(a,t 2 \a,ti) against another uniformly distributed random number R\ between [0, 1]. 

If P(a,t 2 \b,t\) or P(a,t 2 \a,ti) > Ri then the value of the noise at t 2 is x n+2 = a else 
we accept x n+2 = b. By repeating the procedure we can generate a sequence of random 
numbers rj(t) switching between two values a and b. 

It is important to note that the time interval between the two steps is always fixed and 
is equal to At which is much smaller than the correlation time (r) of the noise (At <C r). 

Fig.l(a-b) shows two illustrative dichotomous noise profiles for several values of r. A 
close look into the figures suggests that with increase of correlation time, residence time 
of a particular state increases on an average. In order to check the numerical accuracy 
of the method of generation of dichotomous noise rj(t) we first calculate its first moment 
(r](t)) satisfying the condition (|2.17j) . The typical ensemble averaging is carried out over 
long time of 100 x 5000 time steps as well as over 1000 sample paths. The ensemble average 
(r](t)) is close to zero (less than ~ 10 -6 as determined numerically). In Fig.2 we calculate 

the normalized auto-correlation function ^ w^a? - of the noise for a = 6.0, b = —4.0 for 

{nit) 2 ) 

several given values of r(0.5, 2.0 and 5.0). The values of r are compared to the correlation 
time obtained by the corresponding numerical fitting curves. The agreement is found to be 
satisfactory. 

III. RESULTS: (A) APPLICATION TO NOISE INDUCED CURRENT 

Having explored the method of generation of dichotomous noise we now proceed to solve 
the Langevin equation (|2.1j) or its discretized version (|2.fij) and (|2.7jl . The typical sample 
paths for the particle for several values of noise strength Q are plotted in Fig. 3. It is 
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apparent that while for symmetric noise (9 = 0) the system does not feel any additional 
load, asymmetry in the dichotomous noise gives rise to an external tilt depending on the 
strength of noise. In what follows we note that this asymmetry induced tilt in the overall 
potential is particularly important for generation of current when the potential is symmetric. 
The important relevant quantity for our present study is the steady state average velocity 
which can be defined as 



T is the total time over which the displacements are averaged over for a particular sample 
path. The second averaging is over the number of sample paths, S. Ensemble averaging 
requires both T and S to be sufficiently large. 

We now consider a cosine potential, i.e., V(x) = cos(x), a symmetric periodic potential 
with periodicity 2n. For the entire simulation work we have chosen the time step of inte- 
gration At = 0.01, number of time step is typically of the order of 10 6 and the number of 
sample paths is around 1000. We have calculated the steady state average velocity as a 
function of correlation time r of the dichotomous noise keeping the asymmetric parameter 
9 = and shown (v) vs. r profile in Fig.4. As expected we observe that the average velocity 
is zero, since the potential as well as the external dichotomous noise are symmetric and 
there can not be any directed motion of the particle. Fig.4 thus serves as a thermodynamic 
consistency check for our numerical calculation. 

In order to exhibit Brownian ratchet effect we now switch on the asymmetric parameter 
9 to a non-zero value so that 9^0. In Fig. 5 we display the variation of current (v) as a 
function of correlation time r for three different values of 9(= —1, —2, —4). As the magnitude 
of 9 increases the peak value of the current increases significantly. When r is small (v) tends 
to be vanishing. This is because in order to have a net directed motion it is required that the 
system must not relax to its equilibrium state instantaneously. On the other hand when r is 
very large, the overdamped system gets too much time to reach its equilibrium so that on an 
average net flows in the left and right directions tend to equalize making a net drift current 
vanishingly small. The system behaves resonantly at some optimum value of correlation 
time t of the dichotomous noise, for which the current is maximum. 

In Fig.6(a-b) we exhibit the variation of current ((f)) as a function of correlation time 




(3-1) 
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t of the dichotomous noise for two different values of asymmetry parameter 9 given by 
solid lines using Langevin dynamics. The results are compared with those obtained using 
Fokker-Planck-Master equation^ with periodic boundary condition. In Fig.7(a-b) we make 
a comparison between (v) vs. Q profile (for two different values of r) calculated by our 
Langevin method and the Fokker-Planck method. The agreement is found to be quite 
satisfactory and provides support to our simulation method based on Langevin dynamics. 
(B) RESONANT ACTIVATION OVER A FLUCTUATING BARRIER 
Our second application of the method concerns resonant activation - the well-known 
resonance effec t 39l40l4:Ll42143144145146 which can be observed in the variation of mean first passage 
time as a function of flipping rate of fluctuation of the barrier height of a double well 
potential. For the present purpose we consider a situation where the fluctuation is dichotomic 
in nature, i.e., the barrier height fluctuates between two values. The fluctuation in potential 
has been found to be important in several kinetic model o 43 ' 44 for chemical reactions. We 
consider the potential of the form 

V(x,t) = U(x) + -x 2 r](t) (3.2) 

where U(x) is a bistable potential (— ^x 2 + |i 4 , A and B being constants) with a barrier 
at metastable point x = and two stable minima of the potential are at x = ±y^. The 
fluctuations in the barrier height of the potential are due to dichotomous noise, rj{t) = {a, b}, 
i.e., the barrier height of the potential fluctuates between two states randomly. Under 
overdamped condition the reaction coordinate x(t) is governed by the following Langevin 
equation 

jx = Ax - Bx 3 + xrj(t) + £(t) (3.3) 

£(£) is the usual thermal internal noise as defined in Eq. (j2.2|) and Eq. (j2.3J) An important 
point of departure from the earlier equation (j2.1)l is the presence of a multiplicative noise 
(in the third term on the right hand side of Eq. (j3.3)l ) which is Gaussian distributed and 
exponentially correlated. Since, the stochastic integrals are different in Statonovich and Ito 
definitions, the drift term corresponding to Langevin equation (j3.3J) is Ax — Bx 3 according 
to Ito's prescription, whereas Statonovich prescription provides a drift term of the form 
Ax — Bx 3 + x. The spurious drift component however is missing in ito's definition. It is 
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therefore somewhat convenient to use Ito's scheme for direct simulation of Langevin equation. 
In what follows we adopt this scheme and use of the predictor-corrector steps as (with 
7=1.0); 

zi(*„+i) = x(t n ) + [Ax(t n ) - Bx 3 (t n ) + x(t n )r}(t n )]At + (2k B TAt)V 2 W n (3.4) 

x(t n+1 ) = x(t n ) + ^[{Ax(t n ) - Bx 3 (t n ) + x{t n )rj(t n )} 

+ {A Xl {t n ) - Bx\{t n ) + Xl {t n )r](t n )}}At + (2k B TAt) 1 / 2 W n (3.5) 

where the averages and variances of the white and dichotomous noises are as given earlier 
in the context of Eq. (|2.1j) 

To analyze the essential features of the activated escape following stochastic dynamics 
over a fluctuating potential barrier we numerically simulate the Langevin equation (j3.4H3.5j) 
using the method for dichotomous noise generation as developed in the preceding section 
with a very small time time step (At = 0.01). In our simulation we follow the dynamics of the 
particle starting from the potential minimum at x — — ^J~^oi the left well till it arrives at the 
barrier top at x = where the particle is removed. The first passage time*£ being a statistical 
quantity due to random nature of the dynamics we calculate the statistical average of the 
first passage time over 2000 trajectories. We chose the parameters of the potential as A = 0.5 
and B = 0.1 for the entire set of calculation. We present our simulation results in Fig. 8 
and Fig. 9 for several values of temperatures and variance (^) of dichotomous noise. From 
both the figures it is clear that the mean first passage time ((T)) passes through a minimum 
as one increases the correlation time (r) of dichotomic barrier fluctuation. It is apparent 
that the mean first passage time responds resonantly with the correlated fluctuations of 
the barrier height. As one increases the strength of the internal noise of heat bath the 
escape rate increases (or the mean first passage time decreases) and the resonance behavior 
is manifested at the very high flipping rate as evident from Fig. 8. When the flipping rate of 
the barrier height is very high (- — ► 0) the system effectively feels an average barrier height 
so the mean first passage time is independent of the flipping rate. This is evident form 
the inset of Fig. 8. Fig. 8 and Fig.9 show that the mean first passage time becomes almost 
independent of the correlation time in the limit of high correlation time or slow flipping rate. 
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IV. CONCLUSION 



We have presented an algorithm for numerical simulation for generating dichotomous 
noise and solution of the Langevin equation. As immediate applications of the method we 
have calculated nonequilibrium fluctuation induced current due to asymmetric, exponen- 
tially correlated dichotomous noise in a symmetric periodic potential over a wide range of 
parameter values and resonant activation rate of barrier crossing in a double well potential. 
The method is compared to the Fokker-Planck equation based numerical technique for solv- 
ing stochastic dynamics on a discrete lattice simulated by a Master equation and is found 
to be complementary and well-suited for calculation of ratchet effect. We emphasize that 
while the overwhelming majority of the treatment of various phenomena including barrier 
crossing dynamics, resonant activation, stochastic resonance using dichotomous noise rely 
heavily on Fokker-Planck-Master equation approach we anticipate that the present method 
of solution of Langevin dynamics by direct simulation of dichotomous noise will be useful 
for various purposes in these issues. The simulation method can be easily extended to other 
noise processes involving, for example, three state jump process and for a wide range of 
correlation time and potentials without much difficulty 
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Figure Captions 



Fig.l: Generation of dichotomous noise profile for a = 6, b = —4 (a) r = 0.05 and (b) 
r = 0.4. 

Fig.2: Plots of normalized correlation function vs. i! for a = 6.0, b = —4.0, for three 
different noise profiles with given correlation times r = 0.5, 2.0 and 5.0. Bold circles corre- 
sponds to fitted curves. 

Fig.3: Plot of sample paths of the particle against time for different dichotomous noise 
strength (Q) and asymmetry parameter (8) with D = 0.2 and r — 1.0. 

Fig. 4: Plot of average velocity (v) with correlation time r of dichotomous noise for 
asymmetry parameter 8 = 0, D = 0.02 and Q = 1.0. 

Fig. 5: Plot of (v) vs. r with different values of 8 for constant white noise strength 
D = 0.02 for constant dichotomous noise strength Q = 3.0. 

Fig. 6: (a) Average velocity {(v)) vs. dichotomous noise correlation time (r) profile is 
compared with Fokker-Planck-Master equation method (solid circle) for 8 = —2.0, Q = 3.0 
and D = 0.02. Inset: Same plot but on a shorter time scale, (b) Same as in Fig. 6(a) but for 
8 = -4.0. 

Fig. 7: (a) Average velocity ((v)) vs. dichotomous noise strength (Q) profile is compared 
with Fokker-Planck-Master equation method (solid circle) for r = 0.5, 8 = —2.0 and D = 
0.02. (b) Same as in Fig. 7(a) but for r = 5.0. 

Fig. 8: Mean first passage time ((T)) vs. correlation time of dichotomous noise r for 
D = 0.1 (box); D = 1.5 (circle) and D = 3.0 (triangle). (The inset:mean first passage time 
((T)) vs. correlation time of dichotomous noise r on a log scale for D = 0.1) for o = 2.0 
and b = 1.0. 

Fig. 9: Mean first passage time ((T)) vs. correlation time of dichotomous noise r for 
several values of variances of dichotomous noise [a = 3.0, b = 2.0 (triangle); a = 2.5, b — 1.5 
(circle) and a = 2.0, b = 1.0 (square)] for D = 0.1. 
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